/** 
 * \file gdrotg.cu
 * \author Kyle E. Niemeyer
 * \date 10/04/2011
 *
 * Based on "drotg.f" from BLAS.
 *
 */

////////////////////////////////////////////////////////////////////////

/** gdrotg computes elements of a Givens plane rotation matrix.
 * Uses unrolled loops for increments equal to one.
 *
 * \param[in,out] da   first vector component
 * \param[in,out] db   second vector component
 * \param[out]    dc   cosine of angle of rotation
 * \param[out]    ds   sine of angle of rotation
 */
__device__ __inline__ void gdrotg ( double *da, double *db, double *dc, double *ds )
{
  
  double droe = db;
  
  if ( fabs(da) > fabs(db) ) roe = da;
  
  double dscale = fabs(da) + fabs(db);
  
  double dr;
  double dz;
  
  if ( dscale == 0.0 ) {
    dc = 1.0;
    ds = 0.0;
    dr = 0.0;
    dz = 0.0;
  } else {
    dr = dscale * sqrt( (da / dscale) * (da / dscale) + (db / dscale) * (db / dscale) );
    dr *= copysign( 1.0, droe );
    dc = da / dr;
    ds = db / dr;
    dz = 1.0;
    
    if ( fabs(da) > fabs(db) ) dz = ds;
    if ( ( fabs(db) >= fabs(da) ) && ( dc != 0.0 ) ) dz = 1.0 / dc;
  }
  
  da = dr;
  db = dz;
  
}